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Abstract 

In this article, we present an analytical direct method, based on a Numerov 
three-point scheme, which is sixth order accurate and has a linear execution time 
on the grid dimension, to solve the discrete one-dimensional Poisson equation with 
Dirichlet boundary conditions. Our results should improve numerical codes used 
mainly in self-consistent calculations in solid state physics. 



1 Introduction 

The one-dimensional Poisson equation, 

d 2 <p 

^2 = ~P> a<x< b, (1) 
with Dirichlet boundary conditions, 

<f>(a) = c u 4>(b) = c 2 , (2) 

plays an important role in many branches of science. Particularly, the Poisson equation 
([!]) is essential in self-consistent calculations in solid state physics [fh. In general, 
we have to solve it numerically many times. Therefore, is vital to have the fastest 
and the most accurate numerical scheme to solve it. In this article, we present a very 
efficient direct method, based on a Numerov [^, Q |^] sixth order numerical scheme, 
to solve the Poisson equation ([I]) numerically. Because of its efficiency and simplicity, 
this new method can be used as a canonical numerical scheme to accurately solve the 
one-dimensional Poisson equation. 
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This article is organized as follows. Our numerical scheme is presented in Sec- 
tion H Its linearization, together with a few discussions, are presented in Section ||. 
Our conclusions are presented in Section 0. 



2 The Numerov scheme 

Let 4>i = (f>(xi) represents the solution of (|IJ) at the i-th point, x, ; = a + (;' - \)h, of an 
equally spaced net of step h — (b - a) IN and dimension N + I. Let also (p (k> represents 
the k-th derivative evaluated at the same point x,. Then we can evaluate the solution 
at the nearest neighborhood points of x, using Taylor series [g], 

0, ±I = ftxt ± h) = <Pt ± <pf h + X -<pf h 2 ± 0(h 3 ). (3) 

The basic idea in the Numerov approach is to eliminate the fourth order derivative 
in the expression 

otiAi + a 2 A 2 = arfi - (y a x + a 2 )p< + (^-Q'i + ya 2 )^ 4) + 0{h 6 ), (4) 



where 

l h 1 k 4 



M = ^{<Pm + 4>i-i) = 4>> - '-jPi + 'j^4>f ] + ■■■ , (5) 

A 2 = \{^ + = + Pi-i) = -Pi + jtff + ■ ■ • , (6) 

to obtain the sixth order three-point numerical scheme 

/i 2 

4> i± i = 2cf>i - cf> m - + 10p ' + ^ 

where we chose ci = 1 and, consequently, a 2 = — /i 2 / 12. In a similar way, we can 
eliminate the third order derivative from 

pyB x + B 2 B 2 = hBtff + (jfr + hB 2 )pf ] + 0(h 5 ), (8) 

where 

Bi = l -{<p M - £ -i) = Hf + j<t>f> + • • • , (9) 

B 2 = - <pf x ) = -\{p M - Pi -,) = -Vf » + • ■ ■ , (10) 

to obtain the fifth order three-point numerical scheme 

1 h 2 
<t> f = — - W ) + -(p, + 1 -p w ), (11) 
1 2/z 6 

for the first derivative of where we chose /?i = 1 and, consequently, B 2 = -h 2 /6. 
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So far, the three-point numerical scheme (|TT]) is an iterative method, i.e., given two 
informations, 0,_i and 0;, we can calculate 0,+]. One difficulty of this iterative method 
is related with the Dirichlet boundary conditions they are known only at end- 
points x\ and xn+\- Thus, we can not initiate our iterative scheme (|TT|). Fortunately, 
the recurrence relation in ( pi] ) is linear with constant coefficients. These two features 
imply we can find an unique solution to it, 

h 2 ' 

0; = (i - 1)0 2 - (i - 2)0! - — Y(i + 1 - j)(pj + 10pj-l + Pj-2), (12) 

7=3 

where 0i = c\ and 02 must be expressed in terms of <Pn+\ = c-i (the Dirichlet boundary 
conditions), 

1 1 \r N+l 

02 = -0/v + i + (1 - -)0i + Xi (A? + 2 _ ^ + 10p -'- 1 + ^- 2) ' (13) 

7=3 

Now we have an analytical sixth order numerical scheme to solve accurately the Pois- 
son equation ([!]) with the Dirichlet boundary conditions (||). 

It should be mentioned that the analytical third order numerical scheme presented 
by Hu and O'Connell [[|, making use of tridiagonal matrices, can also be derived by 
the present approach restricted to the third order, 

i 

0,- = (i - 1)02 - (i - 2)0! - h 1 Yj,i + 1 - j)pj-u (14) 

7=3 

where 

1 1 N+1 

02 = -<Pn + i + (1 - -)0i + h 2 YjN + 2- j) Phl . (15) 

7=3 



3 Discussions 

Although we have found a very accurate analytical direct method to solve the one- 
dimensional Poisson equation with Dirichlet boundary conditions, namely, the sixth 
order Numerov scheme (JlTJ), it has one undesirable feature: its execution time is pro- 
portional to the square of the grid dimension. Fortunately it can be linearized. First, 
we create a vector U, whose components are the partial sums Uj = pi + 10p,_i + p,_2 
{U\ — U2 = 0). Next, we create a second vector V with V, = Vj_i + U, and V\ — V2 = 0. 
We also need a third vector Y with T, = iUj and a fourth vector Z with the complete 
sums Z, = Zj_i + Yj. Using these new vectors, our sixth order Numerov scheme dTT| ) 
can be rewritten as follows, 

h 2 

<Pi = a - 1)02 - q - 2)0i - — [a + m - zj. (i6) 

This numerical scheme has now a linear execution time proportional to five times the 
grid dimension N+l. 
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Let us use a Gaussian density, 



p(x) = e 



-x 2 /4 



(17) 



to verify the accuracy and the efficiency of the non-linear numerical scheme (p"2|), as 
well as the linear numerical scheme (16). The solution for the Poisson equation ([!]), 
along with the boundary conditions <p(-10) = 0i = 1 and <p(+lQ) = (/>n+\ = 2, is 



x i 3 

fax) = — - y/7rxerf(x/2) - 2e _J ~ /4 + - + 10 V7rerf(5) + 2e~ 25 , 



where erf(x) is the error function, 

erf( 



x) = — [\' 2 
Jo 



dt. 



(18) 



(19) 



Figure Ha, shows the execution time as a function of the grid dimension iV + 1 for 
three cases. In one case (the dotted line), the numerical solution was computed by the 
non-linear third order numerical scheme ([l4|). In the second case (the dashed line), the 
numerical solution was computed by the non-linear sixth order numerical scheme (|i"2|). 
In the last case (the solid line), the numerical solution was computed by the linear sixth 
order numerical scheme (|l6|). AtN- 1000, the execution time of the non-linear third 
(sixth) order numerical scheme is approximately 145 (51) times the execution time 
of the linear sixth order numerical scheme. Clearly, we can see that the linearization 
process described above plays an essential role in the present Numerov scheme. 

In order to measure the accuracy of the present Numerov scheme, we can compute 
the Euclidean norm 



\W N \\ = 



N+l 



(20) 



where cD w s tands for the exact solution (|l8]) and O*"' stands for the numerical solution. 
Figure 1(b) shows (right vertical axis) a comparasion between two Euclidean norms 
(j20]): one (dashed line) using the third-order numerical scheme dl4| ) and the other (solid 
line) using the sixth-order numerical scheme (16). Note that, at N = 400, the exact 
Euclidean norm of the third-order scheme is approximately four orders of magnitude 
above the exact Euclidean norm of the sixth-order scheme. Naturally, we can see that 
the sixth-order numerical scheme ( |l6| ) is much more accurate and efficient than the 
third-order numerical scheme (|T2|). Of course, we don't know the exact solution in 
practical applications. In that case, the best we can do is to compute the mean Euclidean 
norm of the numerical solution <D W , 



1 /V+1 / \2 
■ /=1 



(21) 



Thi s me an Euclidean norm can be used as a convergency criterion, as shown in Fig- 
ure 1(b) (left vertical axis). 
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Figure 1: Execution times (a) and exact (||VVjvll) and mean (||<J>jv||) Euclidean norms 
(b) as functions of the grid dimension N + 1 . The exact solution is given in ( Jl8| ) and 
corresponds to the Gaussian density ( |l7| ) with boundary conditions 0i = 1 and (p^+i = 
2 

4 Conclusions 

We have applied the Numerov method to derive a sixth-order numerical scheme to 
solve the one-dimensional Poisson equation (|I]) with Dirichlet boundary conditions. 
The resulting recurrence relations were exactly solved and the corresponding execution 
time was linearized [see (|l6])] in such way to avoid the handling of a dense matrix. 
Therefore, the numerical scheme (|l6]) is both accurate and efficient as illustrated in 
Figure [jj. Moreover, it is extremely ease to implement in any numerical or algebraic 
computer language. As pointed by J. M. Blatt |^], the Numerov method is both a 
three-point method, which implies it is stable, and of highest order, which implies it is 
accurate. All these features make the numerical scheme (|l6|) the canonical method of 
choice for the integration of the Poisson equation ([!]). 
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